-
Notifications
You must be signed in to change notification settings - Fork 1.1k
Incidence Angle Modifier zero instead of np.nan (#338) #339
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
Incidence Angle Modifier zero instead of np.nan (#338) #339
Conversation
Could we have the identical code in each of the three functions to do the same thing? That is, to deal with the negative values and the larger than 90 angles? I don't know whether np.maximum, np.where or np.clip is better, but I would lean toward clip. The SAPM model doesn't like negative angles, the others are fine because of the cosine. That should be handled consistently as well. |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Thanks @jkfm. Can you also make sure that the tests assert than the functions preserve nan values? I think they currently do, but this will prevent someone from inadvertently adding e.g. fillna
to the functions in the future.
@@ -811,6 +812,10 @@ def physicaliam(aoi, n=1.526, K=4., L=0.002): | |||
spa | |||
ashraeiam | |||
''' | |||
zeroang = 1e-06 | |||
|
|||
aoi = np.where(aoi == 0, zeroang, aoi) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Why did you add this line? I understand what it does, but not why we need it now when we didn't previously.
It's not clear to me that we actually need zeroang at all. I think its only purpose is to avoid infs, but maybe those be more cleanly handled with np.where or similar.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
As I understood it, if the input angle is 0, tau will be nan, which is why we substitute it with zeroang.
As iam is defined as tau(theta) / tau(0°), we need zeroang anyway to calculate tau(0°). That's why I moved the variable definition to the beginning of the function.
I guess we could define tau0 as a constant in the function and handle aoi == 0 with np.where after the calculations are done, but I feel like the current function is more intuitive to understand.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
We should keep zeroang and for input theta < zeroang, physicaliam should return tau0. Passing theta=0 through will return nan rather than zero. Line 821 computes tau with the leading factor exp(-KL/cos(theta)) which is why we should return tau0 for small theta.
tau0 is computed awkwardly. I'd replace line 829 - 833 with the limit value which we should have implemented in Matlab:
tau0 = exp(-KL)*(1 - ((1-n)/(1+n))^2)
We used zeroang instead in Matlab which is probably why its in pvlib-python
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
ok, now I see that you've changed the test so that an input of 0 returns 1 instead of np.nan. I should have caught this when I wrote the test.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I am surprised that the model would not work for an angle of exactly zero degrees. Could this have something to do with the equation discrepancy mentioned in the comments further above?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I don't think it has to do with the function discrepancy. The unmodified function for theta_r is not defined for angles > ~40°.
The tau(theta)-function in turn simply has a discontinuity for theta = 0. It seems this was ignored in the paper...
@@ -732,7 +732,7 @@ def ashraeiam(aoi, b=0.05): | |||
|
|||
iam = 1 - b*((1/np.cos(np.radians(aoi)) - 1)) | |||
|
|||
iam = np.where(np.abs(aoi) >= 90, np.nan, iam) | |||
iam = np.where(np.abs(aoi) >= 90, 0, iam) | |||
iam = np.maximum(0, iam) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
is the np.maximum line still necessary?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I think it is, because for angles close to +/-90° the function would return a negative value for the iam (the intercept for b=0.05 is at 87.21°: https://www.desmos.com/calculator/y8wvzrx3b9).
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I agree, the np.maximum line is needed
@adriesse brings up a separate but closely related point regarding the meaning of negative values of aoi in pvlib. I tend to agree that the aoi loss functions should all start with something like |
If there are negative values for AOI, that indicates a problem elsewhere. I’m not in favor of using abs(aoi) to address it.
|
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Thanks @jkfm.
If we're not going to treat negative aoi as being equivalent to positive aoi, then my preference is for nans for both AOI < 0 and AOI > 90. Another option is 0 for either or both of those ranges.
In any case, we should update the Variables and Symbols page with the allowable AOI domain. http://pvlib-python.readthedocs.io/en/latest/variables_style_rules.html
@@ -811,6 +812,10 @@ def physicaliam(aoi, n=1.526, K=4., L=0.002): | |||
spa | |||
ashraeiam | |||
''' | |||
zeroang = 1e-06 | |||
|
|||
aoi = np.where(aoi == 0, zeroang, aoi) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
ok, now I see that you've changed the test so that an input of 0 returns 1 instead of np.nan. I should have caught this when I wrote the test.
I vote for IAM=0 where AOI>=90, which is consistent with its definition and its use in the modeling process.
Where AOI=0 (or very close to 0) in the physical IAM model the formula needs to switch to the limit expression in my earlier comment. That switch should be the use of zeroang.
…Sent from my iPhone
|
I would also prefer IAM to return zero for angles outside -90/90 degrees. |
@@ -732,7 +732,7 @@ def ashraeiam(aoi, b=0.05): | |||
|
|||
iam = 1 - b*((1/np.cos(np.radians(aoi)) - 1)) | |||
|
|||
iam = np.where(np.abs(aoi) >= 90, np.nan, iam) | |||
iam = np.where(np.abs(aoi) >= 90, 0, iam) | |||
iam = np.maximum(0, iam) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I agree, the np.maximum line is needed
@@ -811,6 +812,10 @@ def physicaliam(aoi, n=1.526, K=4., L=0.002): | |||
spa | |||
ashraeiam | |||
''' | |||
zeroang = 1e-06 | |||
|
|||
aoi = np.where(aoi == 0, zeroang, aoi) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
We should keep zeroang and for input theta < zeroang, physicaliam should return tau0. Passing theta=0 through will return nan rather than zero. Line 821 computes tau with the leading factor exp(-KL/cos(theta)) which is why we should return tau0 for small theta.
tau0 is computed awkwardly. I'd replace line 829 - 833 with the limit value which we should have implemented in Matlab:
tau0 = exp(-KL)*(1 - ((1-n)/(1+n))^2)
We used zeroang instead in Matlab which is probably why its in pvlib-python
@@ -831,7 +834,7 @@ def physicaliam(aoi, n=1.526, K=4., L=0.002): | |||
|
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Use
tau0 = (np.exp(- 1.0 * (K*L)) * (1 - (1 - n)^2/(1 + n)^2)
instead of the expression on line 826. I'd switch to tau0 when thetar_deg<zeroing.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I agree with your suggestions, especially after having consulted my D&B book. I would suggest calculating reflectance and absorption separately before combining them in order to clarify the code.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I think I applied both your suggestions correctly. Could you take a look?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
It wasn't quite correct, so I thought it might be easier to communicate it this way:
# angle of refraction
thetar_deg = tools.asind(tools.sind(aoi) / n)
# reflectance and transmittance for normal incidence light
rho_zero = ((1-n) / (1+n)) ** 2
tau_zero = np.exp(-K*L)
# reflectance for parallel and perpendicular polarized light
rho_para = (tools.tand(thetar_deg - aoi) /
tools.tand(thetar_deg + aoi)) ** 2
rho_perp = (tools.sind(thetar_deg - aoi) /
tools.sind(thetar_deg + aoi)) ** 2
# transmittance for non-normal light
tau = np.exp(-K*L / tools.cosd(thetar_deg))
# iam is ratio of non-normal to normal incidence transmitted light
# after deducting the reflected portion of each
iam = ((1 - (rho_para + rho_perp) / 2) / (1 - rho_zero)
* tau / tau_zero )
# angles near zero produce nan, but iam is defined as one
small_angle = 1e-06
iam = np.where(np.abs(aoi) < small_angle, 1.0, iam)
# angles at 90 degrees can produce tiny negative values, which should be zero
# this is a result of calculation precision rather than the physical model
iam = np.where(iam < 0, 0, iam)
# for light coming from behind the plane, none can enter the module
iam = np.where(aoi > 90, 0, iam)
# negative angles are not logical
iam = np.where(aoi < 0, np.nan, iam)
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Thanks for the detailed explanation! So does this refer to reference 2 in the Docstring (I don't have access so cannot check)?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Ref [1] simply describes what is in ref [2].
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I implemented your suggestion except for setting negative angles to np.nan, to keep consistency with the other IAM functions for now...
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Consistency has some appeal...
Another observation: contrary to the statement in ref [1] the default value of 2 mm for the glass thickness is definitely at the thin extreme. The value I currently see most often is 3.2 mm. Previously I believe 4 mm was more common. |
I agree with the default value, but I could not find any sources on it. Do you think, we should change the default nevertheless? |
…' into modelchain_zero_instead_of_npnan
The source for the 2mm default has no sources either... |
If it helps, the model output is essentially insensitive to the thickness parameter. You can double from 0.002 to 0.004 and the IAM changes very slightly. For myself, I'd prefer the default to remain the value in the reference paper, simply for traceability back to literature. I do see merit in changing to 3.2mm, though. |
I suggest keeping the default set to 2 mm. I'll merge tomorrow unless there are objections. |
This PR changes the three IAM-models to return 0 instead of np.nan for aoi >= 90 (ashraeiam, physicaliam) and <0 (sapm_aoi_loss).
I also updated the tests to account for the changes. The test for physicaliam now shows, that angles of -90 and 90 return 0, whereas an angle of 0 returns 1, in accordance to the model description in De Soto et al..
Comments are welcome!